#导入仿真优化需要库
import random
import numpy as np
import matplotlib.pyplot as plt
# 用来正常显示中文标签
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus'] = False
import datetime
import simpy
import pandas as pd



#常熟变量的赋值
pop_size = 500  # 种群数量
max_value = 10  # 基因中允许出现的最大值
chrom_length = 100  # 染色体长度
pc = 0.8  # 交配概率
pm = 0.2  # 变异概率
generation = 300
wood_capacity = 500
initial_wood = 300
dispatch_capacity  = 500
days = 100
hours = 8
total_time = days*hours



#运用simpy对工厂的仿真类
class Guitar_Factory1:
    def __init__(self, env):
        self.wood = simpy.Container(env, capacity=wood_capacity, init=initial_wood)
        self.wood_control = env.process(self.wood_stock_control(env))
        self.dispatch = simpy.Container(env, capacity=dispatch_capacity, init=0)
        self.dispatch_control = env.process(self.dispatch_guitars_control(env))
        self.xuqiu = []
        for i in range(20):
            self.xuqiu.append(np.random.normal(loc=200, scale=10, size=None))
        self.day = 0
        self.store = 0
        self.lastmake = 0
        self.lastwood = 500
        self.woodstoreprice = 1.5
        self.woodcost = 1
        self.dispatchstoreprice = 2.5
        self.dispatchcunchu = 0
        self.woodcunchu = 0
        self.profit = 0
        self.producenum = []




#原材料的控制函数
    def wood_stock_control(self, env):
        yield env.timeout(0)
        while True:
            if env.now % 8 == 0:
                yield env.timeout(8)
                if env.now % 40 == 0:
                    self.woodcunchu = (self.lastwood + self.wood.level) * self.woodstoreprice / 2 * 5
                    self.lastwood  = self.wood.level
                    self.profit-=self.woodcunchu
                if self.wood.level <=122:
                    yield env.timeout(16)
                    yield self.wood.put(175)
                    self.profit-=1000+175*self.woodcost




#生产产品的控制函数
    def dispatch_guitars_control(self, env):
        yield env.timeout(0)
        while True:
            if env.now % 8 == 0:
                self.day+=1
                pro = self.producenum[int(env.now / 8)]
                if (pro > 0):
                    self.store += pro
                    self.wood.get(pro)
                    self.profit -= 500
                    yield env.timeout(8)
                else:
                    yield env.timeout(8)
            if env.now % 40 == 0:
                self.day+=1
                self.dispatchcunchu = (self.lastmake + self.store) * self.dispatchstoreprice / 2 * 5
                self.profit -= self.dispatchcunchu
                if (self.store > self.xuqiu[int(env.now / 40)]):
                    self.profit += self.xuqiu[int(env.now / 40)] * 50
                    self.store -= self.xuqiu[int(env.now / 40)]
                else:
                    if (self.store > 0):
                        self.profit = self.profit + (
                                    self.store * 50 - (self.xuqiu[int(env.now / 40)] - self.store) * 20)
                        self.store = 0
                    else:
                        self.profit -= (self.xuqiu[int(env.now / 40)]) * 20
                pro = self.producenum[int(env.now / 8)]
                if (pro > 0):
                    self.store += pro
                    self.wood.get(pro)
                    self.profit -= 500
                    yield env.timeout(8)
                else:
                    yield env.timeout(8)
                self.lastmake = self.store




#运用对照的仿真类
class Guitar_Factory:
    def __init__(self, env):
        self.wood = simpy.Container(env, capacity=wood_capacity, init=initial_wood)
        self.wood_control = env.process(self.wood_stock_control(env))
        self.dispatch = simpy.Container(env, capacity=dispatch_capacity, init=0)
        self.dispatch_control = env.process(self.dispatch_guitars_control(env))
        self.proschedule = list(np.random.randint(0, 2, size=days))
        self.xuqiu = []
        for i in range(20):
            self.xuqiu.append(200)
        self.day = 0
        self.store = 0
        self.lastmake = 0
        self.lastwood = 500
        self.woodstoreprice = 1.5
        self.woodcost = 1
        self.dispatchstoreprice = 2.5
        self.dispatchcunchu = 0
        self.woodcunchu = 0
        self.profit = 0
        self.producenum = []

    def wood_stock_control(self, env):
        yield env.timeout(0)
        while True:
            if env.now % 8 == 0:
                yield env.timeout(8)
                if env.now % 40 == 0:
                    self.woodcunchu = (self.lastwood + self.wood.level) * self.woodstoreprice / 2 * 5
                    self.lastwood  = self.wood.level
                    self.profit-=self.woodcunchu
                if self.wood.level <=122:
                    yield env.timeout(16)
                    yield self.wood.put(230)
                    self.profit-=1000+230*self.woodcost



    def dispatch_guitars_control(self, env):
        yield env.timeout(0)
        while True:
            if env.now%8==0:
                self.day += 1
                if((self.proschedule[int(env.now/8)]==1)&(self.wood.level>100)):
                    self.store+=100
                    self.wood.get(100)
                    self.profit-=500
                    self.producenum.append(100)
                    yield env.timeout(8)
                else:
                    self.producenum.append(0)
                    yield env.timeout(8)
            if env.now % 40 == 0:
                self.day += 1
                self.dispatchcunchu=(self.lastmake+self.store)*self.dispatchstoreprice/2*5
                self.profit-=self.dispatchcunchu
                if (self.store > self.xuqiu[int(env.now / 40)]):
                    self.profit += self.xuqiu[int(env.now / 40)] * 50
                    self.store-=self.xuqiu[int(env.now / 40)]
                else:
                    if (self.store > 0):
                        self.profit =self.profit+ (self.store * 50 - (self.xuqiu[int(env.now / 40)] - self.store) * 20)
                        self.store = 0
                    else:
                        self.profit -= (self.xuqiu[int(env.now / 40)]) * 20
                if ((self.proschedule[int(env.now / 8)] == 1) & (self.wood.level > 100)):
                    self.store += 100
                    self.wood.get(100)
                    self.profit -= 500
                    self.producenum.append(100)
                    yield env.timeout(8)
                else:
                    self.producenum.append(0)
                    yield env.timeout(8)
                self.lastmake = self.store



    def env_status(self, env):
        yield env.timeout(0)
        global status
        status = pd.DataFrame(
            columns=["datetime", "wood", "store"])
        status[
            ["datetime", "wood", "store"]] = status[
            ["datetime", "wood", "store"]].astype(int)
        while True:
            tmp = {'datetime': env.now,
                   'wood': self.wood.level,
                   'store': self.store,
                   }
            status = status.append([tmp])
            yield env.timeout(1)



#基因编码函数
def geneEncodinganquan(pop_size, chrom_length):
    pop = [[]]
    for i in range(pop_size):
        temp = []
        for j in range(chrom_length):
            temp.append(random.randint(0, 1))
        pop.append(temp)
    return pop[1:]



#计算适应度函数
def countvalueanquan(a):
    for i in range(3):
        env = simpy.Environment()
        guitar_factory = Guitar_Factory(env)
        guitar_factory.proschedule = a
        env.run(until=800)
        pro = guitar_factory.producenum
        env = simpy.Environment()
        guitar_factory = Guitar_Factory1(env)
        guitar_factory.producenum = pro
        env.run(until=800)
        return guitar_factory.profit



#生成种群函数
def calobjValueanquan(pop):
    profitlist  =[]
    for i in range(len(pop)):
        profitlist.append(countvalueanquan(pop[i]))
    return profitlist



#计算种群中个体的适应度函数
def calfitValue(obj_value):
    fit_value = []
    c_min = 0
    for i in range(len(obj_value)):
        if(obj_value[i] + c_min > 0):
            temp = c_min + obj_value[i]
        else:
            temp = 0.0
        fit_value.append(temp)
    return fit_value



#求和函数
def sum(fit_value):
    total = 0
    for i in range(len(fit_value)):
        total += fit_value[i]
    return total



#求和函数
def cumsum(fit_value):
    for i in range(len(fit_value) - 2, -1, -1):
        t = 0
        j = 0
        while (j <= i):
            t += fit_value[j]
            j += 1
        fit_value[i] = t
        fit_value[len(fit_value) - 1] = 1



#获取最优个体函数
def best(pop,fit_value):
    best_fit = max(fit_value)
    index = fit_value.index(best_fit)
    return pop[index],best_fit



#获取最劣个体函数
def worst(pop,fit_value):
    best_fit = min(fit_value)
    index = fit_value.index(best_fit)
    return pop[index],best_fit



#选择函数
def selection(pop, fit_value):
    newfit_value = []
    # 适应度总和
    total_fit = sum(fit_value)
    for i in range(len(fit_value)):
        newfit_value.append(fit_value[i] / total_fit)
        # 计算累计概率
    cumsum(newfit_value)
    ms = []
    pop_len = len(pop)
    for i in range(pop_len):
        ms.append(random.random())
    ms.sort()
    fitin = 0
    newin = 0
    newpop = pop
    # 转轮盘选择法
    while newin < pop_len:
        if (ms[newin] < newfit_value[fitin]):
            newpop[newin] = pop[fitin]
            newin = newin + 1
        else:
            fitin = fitin + 1
    pop = newpop



#交叉函数
def crossover(pop, pc):
    pop_len = len(pop)
    for i in range(pop_len - 1):
        if (random.random() < pc):
            cpoint = random.randint(0, len(pop[0]))
            temp1 = []
            temp2 = []
            temp1.extend(pop[i][0:cpoint])
            temp1.extend(pop[i + 1][cpoint:len(pop[i])])
            temp2.extend(pop[i + 1][0:cpoint])
            temp2.extend(pop[i][cpoint:len(pop[i])])
            pop[i] = temp1
            pop[i + 1] = temp2



#变异函数
def mutation(pop, pm):
    px = len(pop)
    py = len(pop[0])
    for i in range(px):
        for i in range(5):
            if (random.random() < pm):
                mpoint = random.randint(0, py - 1)
                if (pop[i][mpoint] == 1):
                    pop[i][mpoint] = 0
                else:
                    pop[i][mpoint] = 1




def mutation1(pop1):
    px = len(pop1)
    py = len(pop1[0])
    obj_value = calobjValueanquan(pop1)        # 个体评价
    fit_value = calfitValue(obj_value)
    best_individual, best_fit = best(pop1, fit_value)
    pop1.remove(best_individual)
    for i in range(px):
        for i in range(5):
            if (random.random() < 0.5):
                mpoint = random.randint(0, py - 1)
                if (pop1[i][mpoint] == 1):
                    pop1[i][mpoint] = 0
                else:
                    pop1[i][mpoint] = 1
    pop1.append(best_individual)



#基于适应度改进的遗传算法操作（并没有运用）
def adanextgeneration(pop,fit_value,num):
    for i in range(0,int(pop_size-pop_size/num),int(pop_size/num)):
        generation1 = pop[i:int(i+pop_size/num)]
        generation2 = pop[int(i+pop_size/num):int(i+pop_size/num*2)]
        fit_value1 = fit_value[i:int(i+pop_size/num)]
        best_individual1, best_fit1 = best(generation1, fit_value1)
        worst_individual1, worst_fit1 = worst(generation1, fit_value1)
        fit_value2 = fit_value[int(i+pop_size/num):int(i+pop_size/num*2)]
        best_individual2, best_fit2 = best(generation2, fit_value2)
        worst_individual2, worst_fit2 = worst(generation2, fit_value2)
        generation1[fit_value1.index(best_fit1)] = best_individual2
        generation2[fit_value2.index(best_fit2)] = best_individual1
        cpoint = random.randint(0, len(pop[0]))
        temp1 = []
        temp2 = []
        temp1.extend(worst_individual1[0:cpoint])
        temp1.extend(worst_individual2[cpoint:len(pop[i])])
        temp2.extend(worst_individual2[0:cpoint])
        temp2.extend(worst_individual1[cpoint:len(pop[i])])
        generation1[fit_value1.index(worst_fit1)] = temp1
        generation2[fit_value2.index(worst_fit2)] = temp2
        pop[i:int(i + pop_size / num)] = generation1
        pop[int(i +pop_size/ num):int(i + pop_size/ num * 2)] = generation2



#仿真过程代码
results = [[]]  # 存储每一代的最优解，N个二元组
fit_value = []  # 个体适应度
fit_mean = []
starttime = datetime.datetime.now()
pop = geneEncodinganquan(pop_size, chrom_length)
num  = 0
for i in range(generation):
    num+=1
    print("第{a}次遗传算法计算:".format(a=num))
    obj_value = calobjValueanquan(pop)        # 个体评价
    fit_value = calfitValue(obj_value)
    best_individual, best_fit = best(pop, fit_value)
    results.append([best_fit, best_individual])
    bestindividual = [[]]
    selection(pop, fit_value)
    for i in range(0,500,50):
        popnew = pop[i:i+50]
        fit_value1 = fit_value[i:i+50]
        best_individual1, best_fit1 = best(popnew, fit_value1)
        bestindividual.append(best_individual1)
    bestindividual = bestindividual[1:]
    print(bestindividual)
    for i in bestindividual:
        pop.remove(i)
    crossover(pop, pc)		# 交配
    mutation(pop, pm)
    mutation1(bestindividual)
    for j in bestindividual:
        p = random.randint(0,500)
        print(p)
        print(j)
        pop.insert(p, j)
    pop[random.randint(0,400)] = best_individual
    print("遗传算法计算结果:", [best_fit, best_individual])



#绘制结果图像
results = results[1:]
X = []
Y = []
for i in range(generation):
    X.append(i)
    t = results[i][0]
    Y.append(t)
plt.figure(figsize=[20,10],dpi = 100)
plt.plot(X, Y,label = "gene")
plt.title("运用遗传算法的迭代结果")
plt.legend()
plt.show()


